function lle = LEs(f, jac, X0, dt, T, warmup, orth_interval)
    if nargin < 6 || isempty(warmup)
        warmup = T/2;
    end
    if nargin < 7 || isempty(orth_interval)
        orth_interval = 20;
    end

    % 预计算步数
    Ntot   = round(T/dt);
    Nwarm  = round(warmup/dt);
    n      = numel(X0);

    % 初始化
    X = X0;
    % 初始化整体变分矩阵 Φ（n×n），只需第一个方向的增长，后面 QR 会正交
    Phi = eye(n);
    sum_log = 0;
    t = 0;

    for k = 1:Ntot
        % —— 1. RK4 积分主系统 —— 
        k1 = f(t,       X);
        k2 = f(t+dt/2,  X + dt/2*k1);
        k3 = f(t+dt/2,  X + dt/2*k2);
        k4 = f(t+dt,    X + dt*k3);
        X  = X + dt/6*(k1 + 2*k2 + 2*k3 + k4);

        % —— 2. RK4 同步积分变分方程：dΦ/dt = J(X)*Φ —— 
        J = jac(X);  % n×n
        L1 = J * Phi;
        L2 = jac(X+dt/2*k1) * (Phi + dt/2 * L1);
        L3 = jac(X+dt/2*k2) * (Phi + dt/2 * L2);
        L4 = jac(X+dt   *k3) * (Phi + dt   * L3);
        Phi = Phi + dt/6*(L1 + 2*L2 + 2*L3 + L4);

        t = t + dt;

        % —— 3. 正交化与累积 —— 
        if k > Nwarm && mod(k,orth_interval)==0
            [Q,R] = qr(Phi,0);
            Phi = Q;
            % 只累积第一个方向的伸缩
            sum_log = sum_log + log(abs(R(1,1)));
        end
    end

    % 用从 Nwarm 到 Ntot 的真实时间
    eff_time = (Ntot - Nwarm) * dt;
    lle = sum_log / eff_time;
end
